In [1]:
import xarray as xr
import numpy as np
import geopandas as gpd
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings('ignore')
In [2]:
plt.style.use('ggplot')
In [3]:
gpm = xr.open_mfdataset('../2005/*')
In [4]:
imd = xr.open_dataset("_Clim_Pred_LRF_New_RF25_IMD0p252005.nc")

Creating a function to correct timestamps of GPM data and to mask the data w.r.t IMD data¶

Also the coordinates of IMD data has been renamed to short and lower form, i.e., TIME --> time, LONGITUDE --> lon, and LATITUDE --> lat.¶
In [5]:
def masked(gpm,imd):
    '''
    gpm: data array to be masked
    imd: maksed dataarray
    returns: masked gpm dataarray
    '''
    gpm1 = gpm.copy()
    gpm1 = gpm1.transpose("time","lat","lon")
    
    imd1 = imd.copy()
    imd1 = imd1.rename({'TIME':'time',"LATITUDE":'lat',"LONGITUDE":'lon'})
    
    gpm1 = gpm1.sel(lat = slice(imd1.lat.min(),imd1.lat.max()),
            lon = slice(imd1.lon.min(),imd1.lon.max()))
    stdTime = gpm1.indexes['time'].to_datetimeindex()
    gpm1['time'] = stdTime
    gpi = gpm1.interp_like(imd1,method='nearest')
    gp = np.ma.array(gpi, mask=np.isnan(imd1))
    gpp = xr.DataArray(gp,coords=gpi.coords)
    return gpp,imd1

Calling Function¶

In [6]:
gpm_new,imd_new = masked(gpm.precipitationCal,imd.RAINFALL)
In [7]:
fig = plt.figure(figsize=[10,4])
ax = plt.subplot(121)
gpm_new.mean("time").plot()
ax = plt.subplot(122)
imd_new.mean("time").plot()
Out[7]:
<matplotlib.collections.QuadMesh at 0x2aab87902340>
In [8]:
gpm_new.mean(["lat","lon"]).plot()
imd_new.mean(["lat","lon"]).plot()
Out[8]:
[<matplotlib.lines.Line2D at 0x2aab88931340>]
In [9]:
gpm_new.sel(time=slice('2005-07-23','2005-07-28'),lat=slice(17,20),lon=slice(71,73)).mean(["lon","lat"]).plot()
imd_new.sel(time=slice('2005-07-23','2005-07-28'),lat=slice(17,20),lon=slice(71,73)).mean(["lon","lat"]).plot()
Out[9]:
[<matplotlib.lines.Line2D at 0x2aab88bced60>]
In [10]:
correlation = xr.corr(gpm_new.mean(["lat","lon"]),imd_new.mean(["lat","lon"]))
print(correlation.values)
0.9064223612416591
In [11]:
bias = gpm_new.mean().values-imd_new.mean().values
print(bias)
0.08560181
In [12]:
import seaborn as sns
In [13]:
sns.kdeplot(gpm_new.mean(["lat","lon"]),imd_new.mean(["lat","lon"]))
plt.xlim(-5,20)
plt.ylim(-5,20)
plt.show()
In [14]:
india = gpd.read_file('IND_adm1.shp')
india.plot()
Out[14]:
<AxesSubplot:>
In [15]:
india[india['NAME_1']=='Maharashtra'].plot()
Out[15]:
<AxesSubplot:>
In [16]:
maha = india[india['NAME_1']=='Maharashtra']
In [17]:
from shapely.geometry import mapping
In [18]:
gpm2 = gpm_new.rio.set_spatial_dims(x_dim="lon", y_dim="lat", inplace=False)
gpm_2 = gpm2.rio.write_crs("epsg:4326")
In [19]:
maharain = gpm_2.rio.clip(maha.geometry.apply(mapping),crs=maha.crs)
In [20]:
imd_2 = imd_new.rio.write_crs("epsg:4326")
In [21]:
mahaimd = imd_2.rio.clip(maha.geometry.apply(mapping),maha.crs)
In [22]:
plt.figure(figsize=[12,4])
ax = plt.subplot(121)
maharain.mean("time").plot.contourf(cmap='jet',levels=range(0,10),extend="max")
maha.plot(ax = ax,ec = 'k',fc='none')
ax.set_title('gpm')
ax = plt.subplot(122)
mahaimd.mean("time").plot.contourf(cmap='jet',levels=range(0,10),extend='max')
maha.plot(ax = ax,ec = 'k',fc='none')
ax.set_title('imd')
plt.show()
In [25]:
import hvplot.xarray
In [109]:
mrain = xr.Dataset({"GPM":maharain.mean(["lat","lon"]),
            "IMD":mahaimd.mean(["lat","lon"])}).drop({"spatial_ref"})
mrain
Out[109]:
<xarray.Dataset>
Dimensions:  (time: 365)
Coordinates:
  * time     (time) datetime64[ns] 2005-01-01 2005-01-02 ... 2005-12-31
Data variables:
    GPM      (time) float32 0.03969 0.2683 0.02317 ... 0.0008004 0.0002058
    IMD      (time) float32 0.1592 0.02913 0.01238 0.04302 ... 0.0 0.0 0.0 0.0
xarray.Dataset
    • time: 365
    • time
      (time)
      datetime64[ns]
      2005-01-01 ... 2005-12-31
      array(['2005-01-01T00:00:00.000000000', '2005-01-02T00:00:00.000000000',
             '2005-01-03T00:00:00.000000000', ..., '2005-12-29T00:00:00.000000000',
             '2005-12-30T00:00:00.000000000', '2005-12-31T00:00:00.000000000'],
            dtype='datetime64[ns]')
    • GPM
      (time)
      float32
      0.03969 0.2683 ... 0.0002058
      array([3.96857113e-02, 2.68319666e-01, 2.31699403e-02, 3.39797349e-03,
             2.49873934e-04, 8.06945463e-05, 2.75250548e-03, 5.54809458e-06,
             1.44165228e-04, 2.08356604e-01, 4.74836415e-04, 0.00000000e+00,
             0.00000000e+00, 1.85418758e-03, 1.07821315e-05, 1.51827437e-04,
             3.51676624e-03, 1.84919185e-03, 0.00000000e+00, 1.44906435e-03,
             9.35438275e-03, 7.19697997e-02, 1.62785791e-03, 5.86837344e-03,
             7.26523757e-01, 3.93860912e+00, 5.04091835e+00, 2.05576563e+00,
             5.32132101e+00, 6.93780184e+00, 3.46774483e+00, 4.08427455e-02,
             6.46763874e-06, 2.24073441e-03, 0.00000000e+00, 2.27773143e-03,
             1.02405012e-01, 3.17980833e-02, 1.98718533e-03, 1.28989015e-03,
             2.15047359e-04, 1.73203042e-03, 1.72765297e-03, 3.43737302e-05,
             2.09796024e-04, 2.71280034e-04, 2.67394260e-03, 1.62055177e-04,
             3.73753556e-03, 4.46575490e-04, 9.52609465e-04, 1.08078726e-01,
             1.34562612e-07, 1.48418239e-05, 1.42123390e-05, 5.74728881e-04,
             1.90347303e-02, 1.67119890e-01, 3.66899490e-01, 1.93798598e-02,
             1.99880172e-03, 2.20123166e-03, 1.00419689e-04, 4.82203934e-04,
             1.05444185e-01, 2.19872117e+00, 4.89956379e+00, 2.55028200e+00,
             1.30486622e-01, 2.78962706e-03, 5.42390253e-03, 1.59397721e-04,
             8.26047093e-04, 4.47369777e-02, 1.52137317e-03, 1.17925601e-03,
             1.36591296e-03, 1.03758706e-03, 4.63520642e-04, 2.04801816e-03,
      ...
             1.72124119e+01, 9.97306919e+00, 3.12388992e+00, 7.23327935e-01,
             1.26292743e-03, 7.60461995e-03, 1.37421663e-03, 1.39616737e-02,
             9.06320438e-02, 8.81063282e-01, 3.88552117e+00, 9.40586746e-01,
             3.26929402e+00, 2.22377276e+00, 1.87040064e-02, 1.14067795e-03,
             2.06307089e-03, 1.61691885e-02, 2.78360234e-03, 1.69321224e-02,
             1.91653191e-04, 6.24027988e-03, 2.85785133e-03, 2.10830495e-02,
             2.35517982e-05, 0.00000000e+00, 0.00000000e+00, 1.57270569e-03,
             5.26800868e-04, 1.56494486e-03, 3.35138146e-04, 5.01077680e-04,
             3.48366587e-03, 3.41243646e-03, 8.22336180e-04, 8.92389915e-04,
             1.14712724e-03, 5.09963487e-04, 6.42003200e-04, 1.92153510e-02,
             1.87865105e-02, 5.89531064e-02, 1.03946127e-01, 2.17885431e-02,
             3.20183896e-02, 1.23287030e-01, 2.57310838e-01, 2.20940262e-02,
             4.59418632e-02, 1.78265143e-02, 1.84880337e-03, 3.60705284e-03,
             6.16506219e-01, 9.54478688e-04, 8.44590017e-04, 4.65323174e-05,
             1.03898114e-03, 1.33298628e-04, 1.43586390e-03, 2.86740862e-04,
             9.01946391e-04, 6.28384296e-04, 7.13398797e-04, 3.40361468e-04,
             2.02612457e-04, 4.49655374e-04, 2.86145747e-01, 3.05633415e-02,
             6.19020546e-04, 5.50999423e-04, 2.51140620e-04, 0.00000000e+00,
             1.15447084e-03, 9.99754848e-05, 4.05460160e-04, 8.00414418e-04,
             2.05781660e-04], dtype=float32)
    • IMD
      (time)
      float32
      0.1592 0.02913 0.01238 ... 0.0 0.0
      array([1.59161136e-01, 2.91277077e-02, 1.23826936e-02, 4.30210307e-02,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 4.52052653e-02,
             0.00000000e+00, 5.40016592e-02, 1.85982659e-02, 0.00000000e+00,
             0.00000000e+00, 8.75096244e-04, 1.35054246e-01, 1.71297565e-02,
             8.84690508e-02, 7.84839988e-01, 2.51101351e+00, 4.10606909e+00,
             2.26715922e+00, 8.05115795e+00, 1.21528206e+01, 2.21456456e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             3.59843895e-02, 1.55394793e-01, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             2.18394189e-03, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 4.64394651e-02, 1.73857480e-01, 5.45727089e-03,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 1.90282594e-02, 1.38457644e+00, 2.91343284e+00,
             1.98267257e+00, 4.60273996e-02, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
      ...
             2.12370758e+01, 1.19436483e+01, 6.90260172e+00, 2.08053970e+00,
             1.01898003e+00, 1.61684584e-02, 4.17529419e-03, 2.76898872e-03,
             3.02927755e-03, 2.55861990e-02, 7.46361554e-01, 2.38397670e+00,
             3.31265666e-02, 2.38448544e-03, 1.17514946e-01, 2.78783534e-02,
             0.00000000e+00, 2.42851079e-02, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             1.02094561e-03, 0.00000000e+00, 6.81242032e-04, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 1.99327478e-03,
             1.04525383e-03, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 1.81030020e-01,
             3.73734580e-03, 0.00000000e+00, 1.42545402e-02, 0.00000000e+00,
             0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
             0.00000000e+00], dtype=float32)
In [77]:
mrain.groupby("time.season").sum().hvplot.bar()
Out[77]:
In [61]:
mrain.hvplot.density()
Out[61]:
In [65]:
mrain.hvplot.hist(alpha=0.5,grid=True)
Out[65]:
In [111]:
plt.figure(figsize=(10,4))
mrain["GPM"].groupby("time.month").sum().plot(label="GPM_SUM",)
mrain["GPM"].groupby("time.month").std().plot(label='GPM_STD',)
mrain["GPM"].groupby("time.month").var().plot.(label="GPM_VAR",)
mrain["IMD"].groupby("time.month").sum().plot(label="IMD_SUM",)
mrain["IMD"].groupby("time.month").std().plot(label='IMD_STD',)
mrain["IMD"].groupby("time.month").var().plot(label="IMD_VAR",)
plt.legend()
Out[111]:
<matplotlib.legend.Legend at 0x2aad05977eb0>
In [335]:
from sklearn.metrics import mean_squared_error, mean_absolute_error
In [336]:
print("RMSE", mean_squared_error(mrain.GPM,mrain.IMD))
print("MAE",mean_absolute_error(mrain.GPM,mrain.IMD))
RMSE 28.681728
MAE 2.3130205
In [338]:
 
In [354]:
x, y =  mrain.GPM,mrain.IMD
m, b = np.polyfit(x, y, 1)
plt.figure(figsize=[6,6])
plt.plot(x, y, 'o')
plt.plot(x, m*x + b)
plt.xlim(-2,62)
plt.ylim(-2,62)
plt.xlabel("GPM")
plt.ylabel("IMD")
plt.show()
del x,y
In [393]:
x, y =  np.arange(1,366),mrain.GPM.values
m, b = np.polyfit(x, y, 1)
plt.figure(figsize=[6,6])
ax = plt.axes()
plt.plot(x, y, 'o')
plt.plot(x, m*x + b)
plt.xlabel("Time")
plt.show()
del x,y
x, y =  np.arange(1,366),mrain.IMD.values
m, b = np.polyfit(x, y, 1)
plt.figure(figsize=[6,6])
plt.plot(x, y, 'bo')
plt.plot(x, m*x + b)
# plt.ylabel("IMD")
plt.show()
del x,y
In [ ]:
 
In [ ]:
 
In [ ]:
 
In [ ]:
 
In [255]:
import pyart
In [256]:
grid = pyart.io.read_grid("../radar_mum/output/grid_MUM150615000342.nc")
In [258]:
grid.fields.keys()
Out[258]:
dict_keys(['REF', 'VEL', 'WIDTH', 'ROI'])
In [259]:
xg = grid.to_xarray()
In [278]:
xg.reindex({'z'})
In [321]:
ds = xr.Dataset(
    {
        "REF":(["time", "z", "lat", "lon"], xg.REF.values,),
        "VELH":(["time", "z", "lat", "lon"], xg.VEL.values,),
    },
    coords={
        "z" :(['z'], xg.z.values),
        "lon" :(["lon"], xg.lon.values),
        "lat" :(["lat"], xg.lat.values),
        "time" :(["time"], xg.time.values),
    }
)
In [338]:
ds['REF'][0,4].plot(cmap='pyart_NWSRef',vmin=-10,vmax=60)
Out[338]:
<matplotlib.collections.QuadMesh at 0x2aad24e31430>
In [344]:
ds["REF"][0].sel(lat=19,lon=71.5,method='nearest').plot(y='z')
Out[344]:
[<matplotlib.lines.Line2D at 0x2aad24f93af0>]
In [337]:
#gpm_new.to_netcdf("gpm_india_2005.nc")
In [ ]: